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Abstract 

We investigate the critical slowing down of the topological modes 
using local updating algorithms in lattice 2-d CP^~^ models. We show 
that the topological modes experience a critical slowing down that is 
much more severe than the one of the quasi-Gaussian modes relevant 
to the magnetic susceptibility, which is characterized by r^ag ~ 
with z ~ 2. We argue that this may be a general feature of Monte 
Carlo simulations of lattice theories with non-trivial topological prop- 
erties, such as QCD, as also suggested by recent Monte Carlo simula- 
tions of 4-d SU(A^) lattice gauge theories. 



Monte Carlo simulations of critical phenomena in statistical mechanics 
and of quantum field theories, such as QCD, in the continuum limit are 
hampered by the problem of critical slowing down (CSD) p. The autocor- 
relation time r, which is related to the number of iterations needed to gen- 
erate a new independent configuration, grows with increasing length scale 
^. In simulations of lattice QCD where the upgrading methods are essen- 
tially local, it has been observed that the autocorrelation times of topological 
modes are typically much larger than those of other observables not related 
to topology, such as Wilson loops and their correlators, see for instance Refs. 
0^111 • Recent Monte Carlo simulations jHUH] of the 4-0? SU(A^) lattice gauge 
theories (for = 3, 4, 6) provided evidence of a severe CSD for the topolog- 
ical modes, using a rather standard local overrelaxed upgrading algorithm 
(constructed taking a mixture of overrelaxed microcanonical and heat-bath 
updatings). Indeed, the autocorrelation time tq of the topological charge 
grows very rapidly with the length scale ^ = o""^/^, where a is the string ten- 
sion, showing an apparent exponential behavior tq ~ exp(c^) in the range of 
values of ^ where data are available. This peculiar effect was not observed 
in plaquette-plaquette or Polyakov line correlations, suggesting an approxi- 
mate decoupling between topological modes and non-topological ones, such 
as those determining the confining properties. The issue of the CSD of topo- 
logical modes is particularly important for lattice QCD, because it may pose 
a serious limitation for numerical studies of physical issues related to topo- 
logical properties, such as the mass and the matrix elements of the rj' meson, 
and in general the physics related to the broken U(1)a symmetry. 

The above-mentioned results suggest that the dynamics of the topo- 
logical modes in Monte Carlo simulations is rather different from that of 
quasi- Gaussian modes. CSD of quasi-Gaussian modes for traditional local 
algorithms, such as standard Metropolis or heat bath, is related to an ap- 
proximate random-walk spread of information around the lattice. Thus, the 
corresponding autocorrelation time r is expected to behave as r ~ (an in- 
dependent configuration is obtained when the information travels a distance 
of the order of the correlation length ^, and the information is transmitted 
from a given site/link to the nearest neighbors). This guess is correct for 
Gaussian (free-field) models; in general it is expected that r ~ where z is 
a dynamical critical exponent, and z ^ 2 for quasi-Gaussian modes. ^ On the 

^Optimized overrelaxation procedures may achieve a reduction of z, although the con- 
dition z > 1 holds for local algorithms |H]- 
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other hand, in the presence of relevant topological modes, the random-walk 
picture may fail, and therefore we may have qualitatively different types of 
CSD. These modes could give rise to sizeable free-energy barriers separat- 
ing different regions of the configuration space. The evolution in this space 
would then present a long-time relaxation due to transitions between differ- 
ent topological charge sectors, and the corresponding autocorrelation time 
should behave as rtop ~ expF;,, where Fh is the typical free-energy barrier 
between different topological sectors. However, for this picture to become 
more quantitative, one should understand how the typical free-energy bar- 
riers scale with the correlation length. For example, we may still have a 
power-law behavior if ~ In^, or an exponential behavior if Ff, ~ It is 
worth mentioning that in physical systems, such as random-field Ising sys- 
tems [53 and glass models PU] , the presence of significant free-energy barriers 
in the configuration space causes a very slow dynamics, and an effective sep- 
aration of short-time relaxation within the free-energy basins from long-time 
relaxation related to the transitions between basins. In the case of random- 
field Ising systems the free-energy barrier picture supplemented with scaling 
arguments leads to the prediction that r ~ exp(c^^), where ^ is a universal 
critical exponent 

Motivated by the recent results of Ref . j3] , suggesting an exponential CSD 
for the topological modes in A-d SU(A^) lattice gauge theories, we decided to 
investigate this issue in 2-d CP^~^ models [121 US], where we can study in 
detail the dependence of the autocorrelation time on the length scale ^ as iV 
is varied. Since the 2-d CP^~^ models possess interesting properties expected 
to hold in QCD, such as asymptotic freedom and a non-trivial topological 
structure, they have often been used as a theoretical laboratory. In particular, 
their lattice formulation has been considered to check and develop methods 
to investigate topological properties in asymptotically free models, exploiting 
also large- analytic calculations. See, e.g., Refs. [Il]-[2n]- The CSD of the 
topological modes in lattice CP^~^ models, and in particular the behavior 
of the autocorrelation time of the topological susceptibility, has already been 
discussed in Refs. [TH| HH] , where the hypothesis of a strong CSD was put 
forward on the basis of a few rough estimates of r^p for the CP^ model, and 
the fact that for large N, N = O(IOO) say, it was not possible to correctly 
sample the topological sectors. In this paper we present high-statistics Monte 
Carlo simulations using local updating algorithms, such as Metropolis and 
overrelaxed algorithms, obtaining rather accurate estimates of the topological 
susceptibility and its integrated autocorrelation time. The results provide a 
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definite evidence that, under local updating algorithms, the CSD experienced 
by the topological modes turns out to be much more severe than the CSD of 
the magnetic susceptibility, whose autocorrelation time shows a power-law 



behavior Tmag ~ with 2; ^ 2. 

Two-dimensional CP^~^ models are defined by the action 

N 
9 



S 



d X D^zDfj_z, 



(1) 



where z is an A^-component complex scalar field subject to the constraint 
zz = 1, and the covariant derivative = <9^ + iA^^ is defined in terms of 
the composite field = izd^z. Like QCD, they are asymptotically free and 
present non-trivial topological structures (instantons, anomalies, 9 vacua). 
The large- expansion is performed by keeping the coupling g fixed [T^ ITB] . 
A topological charge density operator can be defined as 



27r 



with the related topological susceptibility 



Xt 



(fx{q{x)q{0)). 



We consider the lattice formulation ^ ESI EHl EE! 



(2) 



(3) 



(4) 



where, beside the complex A^-component vector z satisfying zz = 1, the 
complex variable A^,^ has been introduced, which satisfies Xn,ij,^n,fj, = ^] Sl 
is a tree-order Symanzik-improved lattice action ITH] . The correlation 
function is defined as 



G{x) = (TrP(x)P(O)), 



(5) 



where P = z ® z. One can define the magnetic susceptibility Xm and the 
second-moment correlation length ^ from its small-momentum behavior: 



Xr. 



G(o), e 



G'(O) - G{q^) 



Asui\qj2) G{q„ 



(6) 
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where gm = (Stt / L, 0) is the minimum non-zero momentum on a lattice of 
size L with periodic boundary conditions (see Ref. [HI] for a discussion of this 
estimator of the second-moment correlation length). We consider the geo- 
metrical definition of lattice topological charge proposed in Ref. [21], which 
meets the demands that the topological charge on the lattice have the classi- 
cal correct continuum limit and be an integer for every lattice configuration 
in a finite volume with periodic boundary conditions. As a result both the 
topological charge and its susceptibility do not require lattice renormaliza- 
tions. It is given by ^29| 



Q = 77-1"^ [In ^T{Pn+^,+uPn+^lPn) + In TT{Pn+uPn+f,+uPn)] , ^J' ^ ^, 

(7) 

where the imaginary part of the logarithm is to be taken in (— 7r,7r). As 
shown in Refs. jTHl |^ , this definition is effective for sufficiently large values 
of A^, where unphysical dislocations [32J should not affect the continuum 
limit of its matrix elements. The corresponding topological susceptibility is 
obtained by 

xt = \^m. (8) 

where V is the volume of the lattice. 

The autocorrelation function Co{t) {t is the discrete Monte Carlo time, 
where a time unit is given by a sweep, i.e. an update of all lattice variables) 
of a given quantity O is defined as 

Co(t) = ((O(t)-(O))(O(0)-(O))), (9) 

where the averages are taken at equilibrium. The integrated autocorrelation 
time To associated with O is given by 

Estimates of tq can be obtained by the binning method (see e.g. Ref. [SS! for 
a discussion of this method and its systematic errors), using the estimator 



ro = TTT^, (11) 



where Eq is the naive error calculated without taking into account the au- 
tocorrelations, and E is the error found after binning, i.e. when the error 
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estimate becomes stable with respect to increasing the block size n?,. The 
statistical error Atq is just given by Ato/tq = ^j2/nb, where is the num- 
ber of blocks corresponding to the estimate of E. As discussed in Ref. jSS] 
this procedure leads to a systematic error of 0{To/b), where b is the size of 
the blocks. In our cases the ratio To/b will always be much smaller than 
the statistical error, so we will neglect it. Equation (fTT|) can be easily ex- 
tended to the case where the quantity O is measured every rim sweeps, i.e. 
To = ^mE^ / {2Eq) , which is of course meaningful only if rim ^ ^o- 

We performed Monte Carlo simulations for = 10, 15, 21, for which 
the geometrical definition ((Zj) should be effective to describe the topological 
modes relevant to the continuum limit. We measured the magnetic suscep- 
tibility Xm, the correlation length ^, the topological susceptibility Xt, and 
the integrated autocorrelation times of Xm and Xt, respectively Tmag and rtop. 
We considered two types of updating methods: a standard Metropolis and a 
mixed method containing overrelaxation procedures. A summary of our runs 
is reported in Tabled Finite-size effects in lattice CP^^^ models are rather 
large and peculiar, especially at large A^ We performed our simulations 
for lattice size L sufficiently large to guarantee that finite-size effects were at 
most O(10~^) for ^ and smaller than 1% for Xt, i-e. L/C, ^ 10 for A^ = 10, 
L/C > 13 for N = 15, and L/C <^ 15 for A^ = 21. Each run consisted typ- 
ically of a few million sweeps for the smallest values of /?, increasing up to 
approximately 50 million for the largest /3's. 

Let us first consider the results of the standard Metropolis algorithm (50% 
acceptance, 10 hits per lattice variable). Figure^ shows the results for the 
integrated autocorrelation times of the magnetic and topological susceptibil- 
ities obtained for the CP^ model. The autocorrelation time Tmag of Xm is 
in agreement with the expected power-law behavior, i.e. Tmag = cC,^ with z 
slightly larger than 2 (a fit of all data to Tmag = cC,^ gives z = 2.30(5) with 
X^/d.o.f. ~ 0.9). On the other hand, the autocorrelation time rtop of Xt ap- 
pears to increase much faster. In particular, a power-law behavior with z ^ 2 
can be definitely excluded. The data for the largest available C suggest larger 
values of z, i.e. z > 4. Moreover, an exponential ansatz rtop ~ exp(c^^) turns 
out to be well fitted by all data, with 6 ^ 0.3, as shown in Fig. Q 

However, in order to obtain a more precise characterization of the topo- 
logical CSD, data for larger values of ^ are necessary, and this becomes rather 
expensive when using a standard Metropolis algorithm. A substantial im- 
provement is obtained by using a more effective local updating algorithm, 
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Table 1: Summary of the Monte Carlo data. The Metropolis and the mixed 
overrelaxed upgrading methods are indicated by mt and ov respectively. 
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Figure 1: Integrated autocorrelation time of the magnetic susceptibility and 
the topological susceptibility for the CP^ model and Metropolis updating. 
The lines are the results of the fits described in the text, power-law and 
exponential fits for Tmag and rtop respectively. 

constructed by employing also overrelation procedures. At each site the 
upgrading method was chosen stochastically between an overrelaxed micro- 
canonical (80%), an over-heat bath [33] (16%), and the Metropolis algorithm 
(4%) to ensure ergodicity. Some details on the application of the above algo- 
rithms to lattice CP^~^ models can be found in Ref. [18j. Similar mixtures 
are usually employed to obtain effective local updating algorithms for 4-d 
SU(A^) gauge theories. The above mixed algorithm turns out to be much 
more effective than the standard Metropolis. For example for = 10 and 
P = 0.70 (^ = 3.8), we found Tmag ~ 380 and rtop ~ 367 using the Metropolis 
algorithm, and Tmag ~ 7 and rtop ~ 27 using the above overrelaxed updating. 
In addition, the mixed algorithm requires less computer time by approxi- 
mately a factor 2. This allowed us to obtain reliable estimates of rtop up to 
,^ ^ 10 for = 10 with a reasonable amount of computer time. Actually, 
we do not exclude that a further improvement can be achieved by optimiz- 
ing the mixture, since we did not really perform a detailed study of this 
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Figure 2: Log-log plot of the integrated autocorrelation time rtop versus ^, 
for the CP^ CP^^ and CP^^ models, obtained using the mixed overrelaxed 
algorithm described in the text. 

issue. ^ Moreover, we performed simulations for larger values of N, N = 15 
and = 21, which will be useful to understand the behavior of rtop, and to 
compare Xt with the available large-A^ results. The results are reported in 
Tabled 

Using the above random mixture of algorithms, the quasi- Gaussian modes 
are expected to be still characterized by power-law CSD with z ^ 2. This is 
substantially confirmed by our simulations. The data for the autocorrelation 
time Ttop of Xt are shown in Fig. |21 It is already apparent from the log- log 
plot of Ttop versus ^ that the data do not agree with a simple power law, i.e. 
with Ttop ~ c^^, on the whole range of ^ explored by this work. Moreover, 
even assuming an asymptotic power-law behavior that sets at relatively large 
^, values z ^ 2 can be definitely excluded, but substantially larger z are sug- 
gested by the data for the largest ^. In the case of = 10 the behavior of rtop 
looks qualitatively similar to the one found using the Metropolis algorithm. 

Comparing the data of rtop at different but fixed ^, we note that the 

^Optimization of overrelaxed algorithms is discussed in Ref. [TB). 
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quantity A^~^ log^g ''"top seems to converge to a non-trivial large- limit, the 
approach being roughly 0{1/N), see also Fig. |21 This fact is also suggested 
by the following simple picture. Let us assume that the transition from one 
topological sector to the other happens by tunnelling through a potential bar- 
rier. The resulting autocorrelation time is r ~ exp 5*;, (neglecting entropy), 
where 5*6 is the action of the typical configurations that are at the boundary of 
the different topological sectors. Let us also assume that these configurations 
are instanton-like. Since the instanton action is given by Si = N2'k / g{p) P?] 
where p is the size of the instanton and g{p) the running coupling at scale p, 
we should expect that Inrtop = 0{N). Note that the same arguments apply 
to A-d SU(A^) lattice gauge theories, and the estimates of the topological au- 
tocorrelation times for = 3, 4, 6 011] are indeed consistent with the above 
dependence. 

Proceeding further within this instanton picture, one arrives at a power- 
law behavior with z ^ N . Let us further assume that the size of the relevant 
instanton configurations at the boundary of different topological sectors is 
given by p ~ a (see, e.g., Refs. jHEni)- Then, we should have g{a) ~ 1/13 
and using asymptotic freedom ^ ~ exp(27r/3), thus rtop ~ exp S'/ ~ with 
z ^ N . As is already apparent from Fig. |2l a reasonable fit to a simple power 
law 

Ttop = fe^r""^ (12) 

can be obtained only by discarding several data points at the smallest values 
of ^. As Fig. 121 already shows, the fitted value of c^v tends to increase when 
more and more data at small ^ are discarded. We can tentatively determine 
lowest bounds for the power coefficient cat by using only the largest values of 
^. For A^ = 10, using the data for the last three /3-values, i.e. data for ^ ^ 7, 
we obtain cio = 0.51(2), hiQ = 0.02(1). In this case, we also note that there is 
a hint of stability in the results for cio, because a consistent value is already 
obtained taking data for ^ > 5, with an acceptable x^- the case A^ = 15, 
the last three /5- values give C15 = 0.49(4), 615 = 0.02(2). Finally, for A^ = 21 
the data for the last three /3- values give C21 = 0.41(2), 621 = 0.07(4). We note 
that these results for are rather close, consistently with the expectation 
that cn = 0(1) in the large-A^ limit. Consistent results are obtained by 
considering a more general ansatz such as rtop = cin + b^C,'^'^'^ , where we 
also allow for a constant term. Note that the naive guess obtained simply 
by assuming g{p = a) = 1/(3 would be 2; = A^. In conclusion, assuming a 
power-law CSD, this analysis indicates that z > N/2. 
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On the other hand, Fig. |21 is also suggestive of an exponential behavior, 
which emerges naturally from tunnelling through a free-energy barrier whose 
size scales like Such exponential behavior gives a good description of 
the data in the whole range of ^ explored by this work. Therefore, we also 
consider an exponential ansatz 

^ log^o nop = aN + bNC^'', (13) 

where qn, and 9n are 0(1) in the large- limit, and 9n ~ 1/2. We may 
further simplify this ansatz by assuming that 6^ is independent of N, i.e. 
6]\f = 6. A global fit to the data gives 6 = 0.49(2). This result was obtained 
by discarding a few data points at the smallest ^ for each A^, i.e. taking 
the data for ^ > 3 (/5 > 0.70) at = 10, ^ > 2 at = 15 and = 21 
(corresponding to (3 > 0.56 and (3 > 0.54 respectively), in order to obtain 
an acceptable x^/d-o.f. ~ 0.9. In Fig. El we also show the curve obtained 
by fits with 6 = 1/2, for which, using the same data as in the above global 
fit, we obtained aio = -0.174(3) and b^ = 0.163(2) with xV^-cf. ^ 1.1, 
ai5 = -0.179(3) and 615 = 0.178(2) with xVd.o.f. ^ 0.4, 021 = -0.169(2) 
and 621 = 0.181(2) with xVd-o.f. ^ 1.0. 

In conclusion, the CSD of the topological modes in Monte Carlo simu- 
lation employing local updating algorithms turns out to be much stronger 
than the one experienced by quasi- Gaussian modes. This has been inferred 
by comparing the integrated autocorrelation time of the magnetic and topo- 
logical susceptibilities as a function of ^. Their behavior suggests an effec- 
tive separation of short-time relaxation within the topological sectors from 
long-time relaxation related to the transitions between different topological 
sectors. A heuristic explanation can be devised by assuming the presence of 
significant free-energy barriers in the configuration space between different 
topological sectors, with the system changing topology by tunnelling through 
such barriers. An exponential ansatz, i.e. xtop ~ exp(c^^) with 6 ~ 1/2, pro- 
vides a good effective description of the data in the range of ^ where data 
are available. However, the statistical analysis of the available data for rtop 
does not actually allow us to distinguish between an exponential CSD and 
an asymptotic power- law behavior with z > N/2 setting at relatively large 
C,. Power-law behaviors with smaller exponents can definitely be excluded by 
our analysis. Data for larger ^ and/or a better modellization of the Monte 
Carlo dynamics of the topological modes would be needed to further clarify 
this issue. We argue that the severe CSD experienced by the topological 
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Figure 3: The quantity ^ logio ''"top versus ^ The hues show the results 
of the exponential fits Ttop ~ exp{biy^^) with 9 = 1/2; see text. 

modes under local updating algorithms should be a general feature of Monte 
Carlo simulations of lattice models with non-trivial topological properties, 
since the mechanism behind this phenomenon should be similar. This is also 
supported by recent Monte Carlo simulations of 4-d lattice SU(A^) gauge the- 
ories reported in Refs. |31 E] . Indeed, the estimates of autocorrelation time 
tq of the topological charge ^ recently reported in Ref. j3] (measured using 
the cooling technique) showed a rapid increase with the length scale, and 
an apparent exponential behavior tq ~ exp(c^) in the range of values of ^ 
where data were available, for all = 3, 4, 6. We stress again that the CSD 
of topological modes may represent a serious limitation for simulations of 
lattice QCD, in order to study physical issues determined by the topological 
excitations, such as the physics of the rj' meson. Our results suggest that the 
contribution of the correlation time to the total cost of a simulation could be 
higher than is usually assumed, if one wants to sample the different topolog- 

■^In our simulation of CP^~^ models we also measured the autocorrelation time tq of 
the topological charge ITjl. In all cases we found rg/rtop — 2 (more precisely tq/ nop — 2.3 
and rg/rtop — 2.1 respectively for the Metropolis and overrelaxed simulations). Note that 
a simple Gaussian propagation would give rg/rtop = 2. 
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Figure 4: Xti^ versus ^ for iV = 10, 15, 21. 



ical sectors correctly. In particular, it may worsen the current cost estimates 
of the dynamical fermion simulations for lattice QCD, see e.g. Ref. [HH] , 
where it is usually assumed that the autocorrelation time only contributes a 
factor of ^. 

An interesting question would be whether other, possibly non-local, up- 
dating algorithms may eliminate or at least improve the severe form of CSD of 
topological modes. Cluster algorithms turn out not to be effective in CP^~^ 
models [13110]. Instead, as shown in Ref. multigrid Monte Carlo al- 
gorithms achieve a substantial reduction of the CSD of the quasi-Gaussian 
modes relevant to the magnetic susceptibility. It is however not clear if they 
can also accelerate the decorrelation of the topological modes. Let us mention 
here that algorithms based on the simulated tempering method [H] were also 
tried in Ref. ^|, but apparently without achieving a particular advantage. 

Finally let us discuss the results for the topological susceptibility; data for 
the dimensionless quantity XtC^) reported in Fig. |1[ clearly show a plateau for 
the largest values of ^, where one can extract an estimate of the continuum 
limit of Xti"^- We obtain: 



Xti'^ = 0.0175(3) for N = 10, 



(14) 
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Xt{2 = 0.0113(2) for N = 15, 
Xt^^ = 0.0080(2) for N = 21, 



where, prudently, we have taken the typical error on the data at the plateau 
as estimate of the uncertainty. Compatible but substantially less precise 
results for = 10,21 are reported in Refs. [IHl HHI- The estimates (fT^ 
may be compared with the available results obtained in the framework of the 
large- expansion 

We note that the estimates (jl4|) are slightly larger. Actually they suggest a 
0(1/A^^) contribution given by C2/N'^ with C2 — 0.2. Indeed, by evaluating 
the quantity N'^ixti'^ ~ 2^)' using the estimates (HD), one would obtain 
C2 = 0.16(3) for A^ = 10,^ ca = 0.16(4) for A^ = 15, and C2 = 0.19(9) for 
A^ = 21. However, this apparent discrepancy can be easily accounted for 
by a slow approach to the large-A^ regime, and the apparent stability of 
the 0(1/A^^) correction may be only a chance. This point deserves further 
investigation. 
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